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Many proposals for quantum information processing are subject to detectable loss errors. In this 
paper, we give a detailed account of recent results in which we showed that topological quantum 
memories can simultaneously tolerate both loss errors and computational errors, with a graceful 
tradeoff between the threshold for each. We further discuss a number of subtleties that arise when 
implementing error correction on topological memories. We particularly focus on the role played by 
degeneracy in the matching algorithms, and present a systematic study its effects on thresholds. We 
also discuss some of the implications of degeneracy for estimating phase transition temperatures in 
the random bond Ising model. 



I. INTRODUCTION 

Quantum information is delicate. The state of a qubit 
can be corrupted by environmental noise, dissipation, 
or by imperfectly implemented logic operations on the 
qubit. Amongst the most remarkable achievements of 
quantum information science are the discoveries of quan- 
tum error correction (QEC) [TJ[2] and fault tolerant quan- 
tum computation (FTQC) [3], which demonstrate that 
storage and manipulation of quantum information is pos- 
sible, even in the presence of such imperfections. Typi- 
cally, this protection is achieved by redundantly encoding 
the logical qubits in a larger number of physical qubits. 

A physical error may be classified as a computational 
error, in which the state of the qubits remains within 
the computational basis, or as a loss^ in which a physi- 
cal qubit (e.g. a photon) is lost from the computer. More 
generally, any detectable leakage process taking the qubit 
out of the computational basis can be treated as a loss 
error. Importantly, losses are presumed to be detectable, 
and in many qubit implementations, losses are a signifi- 
cant source of noise. 

Previous work by Dawson et al. [4 established thresh- 
olds against simultaneous loss and computational errors, 
finding that FTQC is possible provided pioss ^ 0.3% and 
Pcomp ^ 0.01%. Various studies, considering each error 
process in isolation, have derived much higher thresholds 
for each: Varnava et al. [5] constructed an encoding with 
a loss tolerance threshold of 50%, whilst Knill [6 and 
Raussendorf et al. [7 have demonstrated computational 
error thresholds at the 1% level for FTQC. In recent work 
considering the restricted problem of storing quantum in- 
formation in surface codes [8l|9l[T0], we showed that these 
thresholds can be achieved in the presence of simultane- 
ous losses and computational errors [TP. 

The properties of surface code stabilisers that enable 
this synthesis of the known loss and computational er- 
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ror thresholds for quantum memories are shared by the 
more elaborate schemes developed by Raussendorf et al. 
[HITS] that are capable of full FTQC. We therefore fully 
expect that the methods we describe here will apply di- 
rectly to FTQC schemes, including losses. 

In parallel, a fruitful line of research has established the 
intriguing connection between classical error correction 
and statistical mechanics models [M]. More recently, sim- 
ilar connections between quantum error correction and 
phase transitions in many body systems have been made. 
Notably, Dennis et al. [10^ pointed out that optimally de- 
coding an error syndrome on a surface code was equiva- 
lent to solving the classical 2-dimensional random-bond 
Ising model (RBIM). This observation implies that the 
error correction threshold coincides with the Nishimori 
point of the RBIM, i.e. that the error correction thresh- 
old is ^ 10.9% dSl [16]. Further, numerical studies at- 
tempted to not only implement computationally efficient 
error correction algorithms, based on Edmonds' mini- 
mum weight matching algorithm, but to simultaneously 
make predictions about the phase diagram of the RBIM. 
On the basis of large Monte Carlo simulations, Wang et 
al. [15 established very precise predictions about the lo- 
cation of the zero temperature phase transition in the 
RBIM, namely that Pco = 10.31% ± 0.01%. 

The purpose of this paper is twofold. Firstly, we expand 
on the discussion in |llj of a number of subtleties in the 
implementation of our loss tolerant protocol. We provide 
a more detailed derivation of the edge weights that were 
introduced in [iT. We then discuss the 'phase-diagram' 
in the parameter space of losses and computational errors 
giving more detailed results of our numerical Monte Carlo 
simulations, as well as a new analytic result for slope of 
the phase boundary at small values of the loss rate. We 
also give an exposition of the interactions between error 
chains and percolated clusters of loss on finite lattices 
leading to a deviation from universal scaling as the loss 
rate approaches the percolation threshold. 

Secondly, we discuss the role of degeneracy in match- 
ing algorithms for error correction and for the RBIM. 
In the context of error correction, there are in fact very 
many matchings that attain the minimum weight, and so 
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in this sense the code is highly degenerate. This degen- 
eracy was recently discussed in p/7|, and was implicit in 
the edge weighting introduced in [11]. We derive expres- 
sions for the matching probability, accounting for path 
degeneracy. We show that the degeneracy behaves anal- 
ogously to an entropy term in a quasi-free energy. On a 
square lattice, we provide a combinatorial expression for 
the shortest-path degeneracy of a given matching, and 
we exploit this to improve the error correction threshold 
for the surface code modestly, from 10.3% [15 to above 
10.6%. 

In the context of the RBIM, this phenomenon man- 
ifests itself in the very large degeneracy of the ground 
state manifold. The degeneracy becomes important when 
computing numerical estimates of the locations of phase 
transitions in this, and related models. In particular, we 
show that it is important to fairly sample from the ground 
state manifold in order to make accurate quantitative 
predictions of phase transition thresholds. We also give 
some evidence to suggest that existing implementations 
of matching algorithms do not do this, so may possibly 
suffer a systematic inaccuracy. 
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FIG. 1: (a) Physical qubits (arrows) reside on the edges of a 
square lattice (dashed). A plaquette operator is a product of 
four Z operators acting on the qubits within the plaquette 
(blue), and a star operator is a product of four X operators 
acting on the qubits within the star (red). The logical Z op- 
erator is a product of Z operators acting on qubits along the 
solid vertical (blue) line, and the logical X operator is a prod- 
uct of X operators acting on qubits along the solid horizontal 
(red) line, (b) In the event of a qubit loss, an equivalent logical 
operator Z can be routed around the loss. 



II. OVERVIEW OF THE PROTOCOL 

For the purposes of analysis, the error model we con- 
sider is local and uncorrelated. Each physical qubit is 
lost with probability ^loss- Losses are presumed to be de- 
tectable: a projector onto the computational basis of a 
given qubit, 11^ = |0)i(0| + is an observable indi- 

cating whether the state of the qubit has leaked out of the 
computational basis, without affecting the computational 
state of the qubit if it has not. The remaining qubits are 
subject to independent bit-flip (X) and phase (Z) errors, 
each with probability Pcomp- Both errors are handled in 
an analogous way in the surface code, so here we con- 
fine our attention to X errors, noting that the thresholds 
for Z errors will be identical. Aside from these errors, we 
assume other quantum operations (e.g. encoding opera- 
tions, decoding operations, and syndrome measurements) 
can be implemented perfectly. 

Kitaev's surface codes are defined by a set of stabilisers 
[l8] acting on a set of physical qubits that reside on the 
edges of a square lattice. The stabiliser group is generated 
by plaquette operators, which are products of Z opera- 
tors acting on qubits within a plaquette, Pp = ^^^p 
and star operators, which are products of X operators 
acting on qubits within a star, Ss = (8)jGs depicted 
in Fig. [l] [9 . Star and plaquette operators all mutually 
commute, and the code space is a simultaneous +1 

eigenstate of all star and plaquette operators. If the LxL 
lattice has periodic boundary conditions (a torus), there 
are 2L^ physical qubits, — 1 independent plaquette op- 
erators and — 1 independent star operators. There are 
therefore two unspecified degrees of freedom, so the code 
is capable of encoding two logical qubits, Qi {i G {1,2}). 
A logical Zi operator corresponds to a product of Z op- 



erators along a homologically nontrivial cycle extending 
across the entire lattice, shown in Fig. [l] Likewise, a log- 
ical Xi operator corresponds to a cycle of X operators 
extending across the entire lattice in the conjugate direc- 
tion. Zi and Xi commute with the stabilisers, but are not 
contained within the stabiliser group. 

Note that we can obtain a new logical operator by 
multiplying the original one by a plaquette stabiliser, 
Zi = PpZi. The actions of Zi and Zi on the code space 
are the same, since = PpZi\C) = ZiPp\C) = Zi\C). 

Thus there are many homologically equivalent cycles 
spanning the lattice with which to measure each logi- 
cal qubit operator, as shown in Fig. [ijb). This redun- 
dancy allows us to obtain the loss threshold for the case 
Pcomp = 0: if only a few physical qubits are lost, it is likely 
that each logical operator can still be measured by find- 
ing a homologically nontrival cycle which avoids all lost 
physical qubits. Thus the encoded quantum information 
can be reliably recovered. 

If pioss is too high, there is likely to be a percolated 
region of losses that spans the entire lattice, in which case 
there are no homologically nontrivial cycles with which 
to measure the logical operators. In the limit of large L, 
the boundary between recoverable and non-recoverable 
errors is a sharp transition corresponding to the bond 
percolation threshold on the 2D square lattice [19]: for 
Pioss < 0-5 error recovery can almost surely be achieved, 
whereas for pioss > 0.5 error recovery almost surely fails. 
Notably, this threshold saturates the fundamental bound 
on Pioss imposed by the no-cloning theorem pQ] . 

The case pioss = and Pcomp > has been well stud- 
ied m [TOl [15]. Briefly, physical bit-flip errors can lead 
to logical bit-flip {Xi) errors but not logical phase er- 
rors, and vice- versa. An error chain^ is a set of lat- 
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FIG. 2: (a) Lattice with two lost physical qubits. The phys- 
ical qubits are represented by the green lines between pla- 
quettes. Representative physical qubits are labelled 1...8, and 
representative plaquettes are labelled A, B, C, D, E. (b) Con- 
nectivity of (super) plaquettes, where each edge represents a 
physical qubit. (c) Degraded lattice showing superedges (thick 
lines). Black dots indicate the locations of Ising spins in the 
equivalent statistical mechanics model, (d) Restored lattice 
with zero weight edges (dashed) and edges between elements 
of superplaquettes with irregular weights (thick lines). Also 
shown are the test lines for calculating the homology class of 
C = E -\- E' on the restored lattice. 



tice edges (i.e. physical qubits) where a bit-flip error has 
occurred. E changes the eigenvalues of the plaquette op- 
erators only at the boundary, dE, of the chain to —1. A 
connected chain of errors leads to a single pair of pla- 
quettes with eigenvalue —1 at the endpoints of the chain. 
Measuring the plaquette operators therefore yields infor- 
mation about the endpoints of such connected chains. If 
the logical Zi operator crosses the error chain an odd 
number of times, then the logical qubit suffers a Xi er- 
ror. These errors may be corrected if the endpoints can be 
matched by a correction chain, E\ such that the closed 
chain C = E -\- E' crosses the logical Zi chain an even 
number of times, i.e. is homologially trivial. The error 
rate below which the correction chain E' may be success- 
fully constructed is closely related to the phase bound- 
ary of the random-bond Ising model (RBIM) [151 [21]. If 
Pcomp < PcO = 0.104^, then in the limit L ^ oo, the most 
probable chain, Cmax = E -\- E'^^^, is almost surely ho- 
mologically trivial and recovery succeeds. If Pcomp > PcO^ 
then in the limit L ^ oo, the chain is homologically trival 
only 25% of the time. 



^ Some studies report pcO = 0.1031 15 whilst others report pcO ■ 
0.1049 El. 



We can think of the above results as endpoints of a 
'boundary of correctability': (pioss^Pcomp) = (0.5,0) and 
(0, 0.104), respectively. In the remainder of this paper, we 
demonstrate that toric codes (and planar codes, by ex- 
tension) are robust against both loss and computational 
errors, with a graceful tradeoff between the two. We first 
describe how losses can be corrected by forming new sta- 
biliser generators, which are aggregations of plaquettes 
or stars, called superplaquettes and superstars, respec- 
tively. The superstar and superplaquette eigenvalues then 
reveal the error syndromes, and a perfect matching algo- 
rithm is used to find an error correction chain E\ We 
illustrate the efficacy of the single round error correction 
protocol (i.e. ignoring fault-tolerance considerations) by 
calculating numerically the boundary of correctability in 
the {pioss, Pcomp) parameter space. 

Consider the lattice shown in Fig. |2|a) which is dam- 
aged by the loss of two physical qubits, marked by the 
red crosses. The loss of qubit 3 affects two plaquette 
stabilisers: Pa = Z1Z2ZSZ4 and Pb = ZsZ^ZqZ^, ren- 
dering them unmeasurable. However, the superplaquette 
Pab = ^A^B = ZxZ^Z^Z^Z^Z^ is independent of the 
qubit at site 3, and so is a stabliser. In the absence of er- 
rors, Pab has an eigenvlaue of +1. An error chain ending 
within the superplaquette AB changes the eigenvalue of 
Pab to —1. It follows that the syndrome associated with a 
superplaquette is determined by the parity of the number 
of error chains that cross its boundary. The fact that su- 
perplaquette operators yield syndrome information with 
which to construct an error correction chain, E' , is the 
basis for our loss-tolerant error-correction scheme. 

In general, given any set of lost qubits, we can form a 
complete set of stabilisers on the damaged lattice in the 
following way: for each lost qubit which participates in 
neighbouring (super) plaquettes Pg and P^, we form the 
superplaquette operator PgP^, which is independent of 
Zq. In the same way, we form superstar operators from 
products of star operators. As discussed earlier, we can 
also form new logical Xi and Zi operators by deforming 
the original logical operators to conform to the bound- 
aries of newly formed superplaquettes. 

We note that in Fig. [2|a), there is a damaged pla- 
quette operator Z^ = Z^Pa = Z1Z2Z4 (or, equiva- 
lent ly ZsPs = Z^ZqZy) associated with the lost qubit 
3, which commutes with all the newly formed stabiliser 
generators on the damaged lattice, but whose eigen- 
value, ±1, is indeterminate. Likewise, there is a dam- 
aged star operator Xj = X^X^Xs with indeterminate 
eigenvalue that also commutes with the new stabilisers 
on the damaged lattice. Having indeterminate eigenval- 
ues, these mutually- ant icommuting, damaged operators 
form a two-dimensional degree of freedom in an unspec- 
ified state. They therefore describe a completely mixed 
junk qubit, j, which is a consequence of the entanglement 
between the lost qubit and the remaining qubits [22 . 
Since Zj and Xj each commute with the new stabilisers, 
and with the deformed logical operators, the junk qubit 
is therefore in a product state with the logical qubits: 



\i^){^\qi ^ I^)(^lg2 and so the loss does not affect 

the logical qubits. 

When analysing the pattern of syndromes on the pla- 
quettes and superplaquettes, we construct a new graph, 
depicted in Fig. |2|b), in which a (super)plaquette is rep- 
resented by a node, and (super)plaquettes share a bond 
on the new graph wherever they share a physical qubit 
in common. Thus Pab and Pqd share the qubits 2 and 5, 
and this is represented as two edges between the super- 
plaquette nodes AB and CD. 

The error correction syndrome dE arising from an er- 
ror chain on the graph in Fig. [2|b) is determined by the 
(super) plaquettes that have an eigenvalue of —1. To cor- 
rect the errors, we follow the procedure described in pre- 
vious work [lOl [15] to find the most likely error chain 
giving rise to dE. 
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FIG. 3: Supposing the chain of errors E to have been pre- 
ciously laid out, the probability of then obtaining the path 
C by further random flips, given by Eq. ([2]), is calculated is 
calculated with reference to this illustration. The four black 
circles represent the regions over which the four products are 
taken in the first line of Eq. ([2| . 



A. Super-edge Weights 

Calculating the probability of a given error chain is 
complicated by the presence of losses. In the case where 
Pioss = the probability of an error on a qubit, I = 
{z, j}, between two neighbouring plaquettes i and j, is 
uniform, = Pcomp- With losses, superplaquettes may 
share common physical qubits, as shown in Fig. [2|b). In 
this case, the superplaquettes {AB, CD} have two qubits 
in common, 2 and 5, each of which might suffer an error. 
A non-trivial syndrome arises only if either qubit 2 or 
qubit 5 suffered an error, but not both. 

By extension, for a pair of neighbouring superplaque- 
ttes, i = {P, P'} sharing n£ physical qubits, a non- trivial 
syndrome arises only if there are an odd number of errors 
on the n£ qubits, which happens with probability 



m odd 

= (1 - (1 - 2p,„^j,)"^)/2. 
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We therefore degrade the graph shown in Fig. |2[b) , re- 
placing multi-edges (i.e. several shared physical qubits) 
whose error probabilities are uniform, with single su- 
peredges whose error rates depend on the number of phys- 
ical qubits shared between neighbouring superplaquettes. 
This degraded lattice is shown in Fig.|2|c), in which there 
are no multi-edges, but the error probabilities are no 
longer constant. Thus, the edge weights account for the 
degeneracy of possible error paths between neighbouring 
superplaquettes. In subsequent sections, we will discuss 
the general issue of degeneracy on an arbitrary degraded 
lattice. 



B. Chain probabilities 

Given a measured error syndrome indicating the 
boundary, dE, of errors on the degraded lattice, the prob- 
ability of a given error chain E' may be computed. For 



completeness, and to define notation, we recapitulate the 
derivation of Wang et al. [15^ (eq 15), generalised to in- 
clude non-uniform edge weights. Given an error chain E, 
we need a probability distribution for correction chains 
E' that satisfy dE' = dE, that is we need P{E'\E). Since 
C = E' -\- E (which happens to be a closed chain) it 
follows that p\e'\E) = P{C\E). The quantity P{C\E) 
should be interpreted as "Given the chain E of flipped 
edges (i.e. qubits), what is the probability that indepen- 
dently flipping more edges results in only the chain C". 
The set of edges that would need to be randomly flipped 
to transform E into C is £ e {C nE)U (C HE). All other 
edges should remain unflipped, as depicted in Fig.jS] The 
probability of this occurring is 
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In the last line, we have used the identity u 
which follows from the fact that E' = C -\- E. Note that 
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FIG. 4: (top) Statistical sampling of pfaii vs. Perr- (bottom) Same data scaled according to x = (p — Peo)iv'/"°, as in Harrington. 
Then pfaii{x) = a + 6x is fitted to the scaled data. Sample standard-error bars are shown. Red circles are for L = 16, green 
squares are for L = 24 and blue diamonds are for L = 32. 



nf = {1 — uf)/2, and that the definition of Ji differs 
slightly from that in ^15^ . 

On this degraded lattice with non-uniform edge 
weights, we may therefore assign a probability for any 
hypothetical chain E' = E ^ C. This probability, which 
is conditioned on the syndrome dE^ is [TOl [T5] 



P{E'\dE)=J\f\[e^'''^^ 



(4) 



Following [15 , satisfies the closed chain constraint 
n^Gs '^F ~ -'- every (super)star 5, which can be 
'solved' by introducing the Ising spins, and CFj on the 
dual lattice such that u^^^ = (TiCjj. Then 



P{E'\dE) = M\[e^^^<''^''^ , 

{ij) 



(5) 



and this probability is generated by the partition function 
for the RBIM 



^ (ij) 



(6) 



at /3 = 1. There are thus two additional sources of ran- 
domness in this model compared to the standard RBIM: 
(a) the lattice connectivity is irregular, and (b) the mag- 
nitude of the Ising couplings, Ji are irregular. However, 
the probability distribution from which uf is selected 
and the magnitude of the coupling are both functions 
of the probability p£. Whilst this varies throughout the 
lattice, each edge satisfies the a local Nishimori condition. 



A ground state of this irregular RBIM (/3 oo) can 
be found by maximising P{E'\dE) over possible paths 
E' . The path E' that maximises P{E'\dE) is the same 
as that which maximises hiP(yE'\dE): 
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For the regular lattice, p£ = p, and so this amounts 
to finding the shortest path. For a lattice with irregu- 
lar weights (as we have), we should weight each edge by 
ln((l — Pi)/p()^ and find the path with the minimum ad- 
ditive weight. 

The most likely error chain is given by 
max^;/ P{E'\dE). This is equivalent to 
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This minimisation may be accomplished using Edmonds' 
minimum- weight, perfect-matching algorithm [23 . For 
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the case pioss = 0, this simply minimises the total 
metropolis length of the matching path, and is the same 
procedure implemented in previous studies [lOl [15]. In 
the case where pioss > 0, the edges do not have uni- 
form weight, since p£ depends on the number of physical 
qubits, n^, shared between adjacent superplaquettes. 



C. Correctable Phase Diagram 

In order to make quantitative estimates of the tolera- 
ble error rates in the presence of loss, we take adopt two 
approaches. Firstly we perform Monte Carlo simulations 
to provide statistical estimates of the boundary for all 
values of pioss- Secondly, we provide a semi-analytic es- 
timate for the threshold value of Pcomp as a function of 
Pioss ' These two approaches are in good agreement. 



1. Monte Carlo Simulations 

For the purposes of simulation, it is easier to deter- 
mine homology classes on a square lattice, rather than 
the degraded lattice, exemplified in Fig. [2|c). We there- 
fore restore the square lattice by dividing superplaquettes 
into their constituent plaquettes in the following way: (1) 
an edge between two plaquettes within a single superpla- 
quette is assigned a weight of zero, (2) an edge between 
plaquettes in two neighbouring superplaquettes is given 
the weight of the superedge in the degraded lattice, as 
illustrated in Fig. [2|d). These transformations do not 
change the minimum weighted-distance between any pair 
of syndromes, and so a minimum- weight perfect matching 
on the restored lattice is also a minimum-weight perfect 
matching on the degraded lattice. Determining the ho- 
mology class is then accomplished by counting crossings 
of vertical and horizontal test lines in the dual lattice. 

In order to test the efficacy of our loss-tolerant error 
correction scheme, we generate random losses on a peri- 
odic lattice with rate pioss- On the remaining edges we 
generate a set of errors, with rate Pcomp- Depending 
on the distribution of losses, we assign edge weights ac- 
cording to Eq. ([T]) to edges in the restored square lattice. 
The syndrome, dE, is calculated, and applying Edmonds' 
algorithm to dE yields the maximum-likelihood error cor- 
rection chain, E' . The homology class of the chain E-\-E' 
then determines whether error correction was successful. 

For a given value of pioss we simulate the protocol 
for different values of Pcomp on different sized lattices, 
L = 16, 24 and 32. For each value of Pcomp and L, the 
failure rate, pjaii^ is calculated by averaging over 10^ tri- 
als. Following [15 , we seek a threshold, p^omp (depending 

as 
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on Pioss), such that if Pcomp < p\ 
L ^ oo. Conversely, if 
as L ^ oo. 

Results of simulations for different values of pioss are 
shown in Fig. [4j For each value of piossi the error rate 
is scaled according to the same prescription as given by 
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FIG. 5: Correctability phase diagram. The shaded region is 
correctable in the limit L ^ oo. Red points are the crossing 
point of Pfail calculated on L = 16 and L = 24 lattices, green 
are the crossing of pfaii (versus Pcomp) on L = 16 and L = 32 
lattices, and blue are the crossing point of Pfaii on L = 24 and 
L = 32 lattices. The curve is a quadratic fit to the points in the 
interval < pioss < 0.4 for which finite size effects are negligi- 
ble. It extrapolates through the point {pioss, Pcomp) = (0.5, 0). 
The dashed line is the linearised phase boundary for pio ss 
much less than the percolation threshold, given in Eq. (10). 



Wang et al.: x = {p — The scaled results 

generally fall on the same linear relationship, Pfaii = a-\- 
bx^ also shown. The parameters p^^J^^, z/q, <^ and b with 
Pioss and are determined by least squares fitting. If there 
is universal scaling behaviour, then p^^mp independent 
of L, however at large values of pioss, universal scaling 
breaks down, and the point at which the Pfaii curves 
cross depends on L. 

Fig. [5] is the central result in this paper, and shows Pco 
as a function of pioss- Since we simulate three different 
values of L, there are in principle three different candi- 
dates for Pcomp 1 indicated by the different coloured/sized 
points. For pioss ^ 0.4, these all coincide, so p^^mp is well 
defined. The failure rate obeys the universal scaling law 
Pfail = f[{Pcomp - Pcomp)^^^''°] with Scaling exponent, 
z^o, in the range 1.4 to 1.5, consistent with the RBIM 
universality class p!5] . 

Fitting a quadratic through the points in the interval 
< Pioss < 0.4 yields a curve that extrapolates through 
the point (0.5,0). This curve represents the boundary of 
correctability: if (pioss^Pcom^p) is in the shaded region then 
error correction almost surely works, as L ^ oo. Impor- 
tantly, this boundary passes through the known bounds 
at Pioss = and 0.5, demonstrating that the protocol is 
very robust against loss. 



2. Linearised approximation for small pioss 

If Pioss is much less than the percolation threshold, so 
that losses are very sparsely distributed, then to a very 
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FIG. 6: Fit parameters, a, b and z/q, for data collapse in Fig. [Ij^bottom). 



good approximation, the only superplaquettes that arise 
are those consisting of two neighbouring plaquettes, and 
contain 6 physical qubits. The fact that these superpla- 
quettes contain more physical qubits than the original 
plaquettes means that they have a higher chance of suf- 
fering a syndrome error. Therefore, as pioss increases, the 
the effective rate of syndrome errors (averaged over pla- 
quettes and superplaquettes) also increases. We use this 
observation, that for a fixed computational error rate, a 
lossy lattice has a higher per-stabiliser error rate than a 
lossless lattice, to calculate the dependence of the error 
correction threshold on the loss rate, for small pioss- 

For small loss rates the probability that a given sta- 
biliser has a syndrome error is, to a good approximation: 

P(stab. error) = P(no loss)P(plaq. error) 

+ P(l loss)P(superplaq. error) (8) 



\o) by differentiating both sides of Eq. ^ 
with respect to pioss- The result is 



a 



Pcomp 



= -2pco(l -Pco)(l - 2pco), 
« -0.148. 

It follows that the linear approximation to the phase 
boundary for small pioss is 



^thr 



Pc 



' Pco -^apio 



(10) 



This relationship is shown as the dashed line in Fig. |5j 

From the quadratic fit to the Monte Carlo results in 
Fig. ^ a is found to be -0.154 ± 0.0033 (la). This is 
in good agreement with the semi-analytic estimate given 
above. 



where 

P(no loss) 
P(plaq. error) 
P(superplaq. error) 



1 - P(l loss) = (1 



■PlossY 



'^Pcomp^^ /2 



are the probabilities that an odd number of physical 
qubits suffer an error on a plaquette (consisting of 4 
qubits) and on a superplaquette (consisting of 6 qubits), 
respectively^. 

The effective single qubit error rate, Peff^ on a lossless 
lattice that duplicates the per-stabiliser error rate of the 
lossy lattice, Eq. ([8|, therefore satisfies 



(1 - (1 - 2peff)f/2 = P(stab. error). 



(9) 



From [TF , the threshold for this effective error rate on the 
lossless lattice is Pco = 0.103. Thus, setting p^ff = Pco in 
the LHS of Eq. ([9| yields an implicit equation relating the 
threshold error correction rate, Pcomp? Pioss- Treating 
Pioss as the independent parameter determining the error 

so Pcomp ~ PcompiPloss) 



threshold, so p^^^Lr. = pi^L,r){pioss)i we compute the slope 



The latter two are specific cases of Eq. 



D. Superplaquette Percolation and Finite Size 
Effects 



For Pioss ^ 0.425, the universal scaling law breaks 
down, and the points in Fig. [5] lie below the quadratic 
extrapolation (though also reaching the point (0.5,0)). 
This is attributed to the fact that for pioss ^ 0.425, the 
largest superplaquette on an L < 32 x 32 lattice occupies 
approximately half of the lattice sites [24], so finite size 
effects dominate. This can be seen in Fig.[6j in which the 
fit parameters a, b and change precipitously for large 
values of pioss 

As Pioss approaches 0.5 (for fixed L), the size of the 
largest superplaquette becomes comparable to the size of 
the lattice, so it is reasonable to expect that the scaling 
behaviour at small pioss breaks down. The size of the 
largest superplaquette, jULiPioss)^ on an L x L lattice is 
closely related to the size of the largest percolated cluster, 
discussed by Bazant [2T. We note that percolated cluster 
of losses may contain 'islands' of qubits that were not 
lost. For our purposes, these isolated islands within a 
cluster of losses neither assist in constructing deformed 
logical qubit operators, nor in creating superplaquettes. 
Therefore, in analysing the size of percolated clusters, we 
treat these islands as if they are lost. 
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FIG. 7: (left) Scaling law for the size of the largest superplaquette on an L x L lattice, /xl is the number of sites in the largest 
superplaquette, Sp is the number of sites within a 'percolated-correlation length', Sp — \p — Pperc\^^ , where Pperc — 1/2 is the 
percolation threshold for the lattice, v — 4/3, is the correlation length scaling exponent, D — 91/48 is the fractal dimension of 
the largest percolated cluster and d — 2 \s the dimensions of the lattice. Points are statistical samples on different lattice sizes 
ranging from L = 2 to 64, and the solid curve is a fitted scaling law, given by $(x) = 0.116 log(7.74x^/^ + 1) (adapted from 
[24 ). (right) Solving iil{p scale) = iv^/2 for p scale as a function of L. When pioss = P scales the largest super-plaquette occupies 
half of the lattice. 



Fig. [rj^left) shows fiL versus L^, with fiL and are 
scaled by Sp = \p—Pperc\^^ and Sp^^ respectively (where 
ly = 4/3 is the percolation length scaling exponent [El 
[25] , D = 91/48 is the fractal dimension of the largest 
percolated cluster [24] and d = 2 is the dimensions of the 
lattice, and for the purpose of this discussion, p is the 
loss rate). The results collapse onto the universal scaling 
law 



Ml 

Sp 



d/D 



(11) 



where ^{x) = xlog(^ x^^^ + 1), shown as the solid line 
with fit parameters x = 0.116 and ^ = 7.74. 

Intuitively, when the largest superplaquette is compa- 
rable to the size of the lattice, the effective diameter of 
the lattice is reduced, leading to an enhanced probability 
of attempting to match syndromes with homologically 
non-trival paths. This 'interaction' between the loss rate 
and the computational error rate leads to non-universal 
scaling, and degrades the performance of the error cor- 
rection. 

Fig. [Tj^right) shows the loss rate at which the largest 
superplaquette occupies half of the lattice, i.e. Pscaie is 
defined such that liLiPscaie) = ^^/2. Whenp^oss > Pscaie, 
the scaling law expected for small values of pioss breaks 
down. When L = 16, Pscaie = 0.42, which coincides with 
the point at which the red data points in Fig. [5] depart 
from simple quadratic behaviour. When L = 24, Pscaie = 
0.44, which coincides with the point at which the blue 
data points depart from the simple quadratic behaviour. 

As L ^ oo, Pscaie ^ 1/2, which is the percolation 
threshold. Therefore, on larger lattices, the anomalous 
region, pioss > Pscaie becomes smaller. 



III. USING DEGENERACY TO IMPROVE THE 
ERROR CORRECTION THRESHOLD 



The role of degeneracy in improving the performance 
of error correction algorithms has recently been noted by 
Duclos-Cianci and Poulin [17]. For a given configuration 
of the error syndromes, there may be a large number 
of degenerate minimum weight matchings. Fig. [8] exem- 
plifies this for a configuration in which the eigenvalues 
of the four numbered plaquette have flipped. There are 
two degenerate matchings of the numbered plaquettes for 
this case: {{1, 2}, {3, 4}} or {{1, 3}, {2, 4}}. The former 
matching has only a single matching path^ shown in green. 
The latter matching has nine degenerate matching paths: 
there are three paths of minimum weight between 1 and 
3, and likewise between 2 and 4. The higher degeneracy 
of the latter matching indicates that it is more likely: 
there are more minimum-weight error chains that could 
have generated it. 

When trying to find the most likely configuration of er- 
rors given a certain error syndrome, this example demon- 
strates that one should account for the fact that some 
matchings have a higher degeneracy than others. We now 
show that the path matching degeneracy can be handled 
by a small modification to algorithms based on Edmonds' 
perfect matching. 

For a given matching, M, the degeneracy of the match- 
ing. Dm, is the number of shortest paths that are consis- 
tent with the matching. In the case of a square lattice this 
degeneracy is simple to calculate, by considering the in- 
dividual pairings within the matching M = {..., m^b, ...}, 
where niab is the pair of syndrome sites {a, 6}. For each 
pairing, rriab there are hm . horizontal edges and Vm u 
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path probability is then 



FIG. 8: There are two minimum weight matchings of the num- 
bered nodes: {{1, 2}, {3, 4}} or {{1, 3}, {2, 4}}. The former 
has only one minimum weight matching path (green lines), 
whilst the latter has nine, (including the light blue line). 



vertical edges in any shortest path between the sites a 
and b. The set of all shortest paths between a and b is 
therefore the number of permutations of a string of hmab 
horizontal steps and Vmab vertical steps 



'^ruab'^rnab' 



The matching degeneracy is then 
Dm= II 



(12) 



(13) 



meM 



For the example shown in Fig. [Sj the matching 
{{1, 3}, {2, 4}} has pairings {1, 3} and {2, 4}. For the pair 
{1,3} there is hi^s = 1 horizontal edge and ^1^3 = 2 
vertical edges on any shortest path between the plaque- 
ttes. It follows that there are ^{1,3} = (2) = ^ differ- 
ent shortest paths between 1 and 3. The same is true 
for the other pair {2,4}: /i2,4 = 1 and 'L'2,4 = 2, so 
^{2,4} = 3. The degeneracy of the matching is then 

^{{1,3}, {2,4}} = ^{1,3}^{2,4} = 9- 

We wish to find the matching that maximises the 
matching probability P(M). There are Dm shortest 
matching paths, {E[^ ...,£^^^}, and since they are all of 
the same length, they each have a probability 



P{M) = ^ P{El), 
= DmP{E[), 



(14) 



Maximising P{M) is therefore equivalent to finding the 
matching that minimises 



G(M)= J2 {dm-lnD„ 



(15) 



meM 



This minimisation may be accomplished using the 
same methods described earlier, using Edmonds' perfect 
matching, where the distance between two sites a and b 
on the lattice is reduced by the amount In Dm^t, to ac- 
count for the number of degenerate paths between a and 
b. 



The form of Eq. (14) suggests an analogy between the 



degeneracy and the entropic contribution to a free energy, 
G, in a partition function. A temperature scale does not 
appear explicitly in this expression: it is implicitly the 
'partition function' at T = 1.^ It follows that the match- 



ing that minimises Eq. (15) is not necessarily a shortest 
path. 

To make the analogy with a free energy more explicit, 
we modify Eq. ( 15 ) to 



Gr{M) 



E 

meM 



r In Dr, 



(16) 



interpreting r as a temperature-like parameter. When 
minimising Gr in the limit r <C 2J ~ 4 (i.e. much less 
than the RBIM grain boundary tension) excitations out 
of the degenerate ground state manifold are suppressed, 
so only microstates within the ground state manifold con- 
tribute to the entropy. 

In any case, we can take r to be a free parameter with 
which to tune the performance of the matching. We can 
think of r as parameterising a class of matching algo- 
rithms. In order to study the properties of this class of 
matching algorithms, we compute the error correction 
threshold, Pco{t) for different values of r. ^ This is shown 
in Fig. [9}^ 

There are a number of conclusions to draw from Fig. [9] 
Firstly, for < r < 1.5, the error threshold is higher than 



= M'e 



= PiE[), 



where dm^i, is the distance between sites a and b. In the 



absence of loss, dr, 



2J{hn 



). The matching 



The analogy is not perfect, since we should include all matching 
paths for T > 0: whilst Eq. explicitly includes all matchings, 
it only accounts for the shortest matching paths consistent with 
a given matching. 

We note that all numerical results we report here in regard to 
the degeneracy are at pioss = 0- 

For each value of Pcomp: the same 10^ samples were used to 
estimate Pfaii for temperatures 0.1 < r < 1.5. As a result, the 
error bars in Fig. [9] are not uncorrelated. Results at r = are 
taken from TS . 
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FIG. 9: The critical error probability as a function of the 
temperature parameter r. The blue points were calculated on 
even lattices with L = 16, 20, 24, 28 and 32. The red points 
were calculated on odd lattices with L = 15, 19, 23, 27 and 
31. The same error instances were used to calculate pc for all 
0.1 < r > 1.5, so the error bars are correlated. The values at 
T = are taken from T5^. The value at r = 2 is shown to 
demonstrate the dramatic drop with increasing r. 



FIG. 10: The fraction of matchings that are minimum dis- 
tance, for different values of r. As r increases, an increas- 
ing proportion of the selected matchings are longer than the 
minimum distance matchings (i.e. ground state matchings), 
but have sufficiently high degeneracy that their free energy 
is lower. All simulations shown were done at Pcomp = 0.1035, 
on different sized lattices, as indicated. Note the disparity be- 
tween the minimum distance matching fraction on even (blue) 
and odd (red) lattices. 



that at r = 0. For small r, the degeneracy in Eq. ( 16 ) car- 
ries a small weight, so that G is maximised with match- 
ings in the ground state manifold: the matching algo- 
rithm chooses a maximally degenerate matching within 
that manifold. This is demonstrated Fig. [lO) where the 
fraction of matchings that are also minimum distance 
matchings is shown for different values of r. For r < 0.5, 
this fraction is essentially unity, so that all matchings are 
also minimum distance matchings. 

The threshold increases with increasing r, and peaks 
around r = 1.4 or so. Interestingly, this is somewhat bet- 
ter than the threshold at r = 1, pertaining to Eq. ([Ts]). 
When r > 1, the matching algorithm no longer sim- 
ply maximises the matching path probability, P(M). 
Instead, it favours more degenerate paths. Finally, the 
threshold decreases as r increases further. This behaviour 
is qualitatively consistent with the proposed re-entrant 
phase diagram for the RBIM, shown in e.g. [15]. 

There is a notable discrepancy between the thresholds 
and minimum distance matching fractions on even and 
odd lattices. For r < 1 even lattices appear to have a 
substantially higher threshold than do odd lattices. This 
discrepancy is much larger than that discussed in [15], 
shown for r = 0. The deviation of the minimum weight 
matching fraction from unity also follows different trends 
depending on the parity of the lattice: on odd lattices it 
begins deviating from unity for r > 0.6, whilst on even 
lattices the deviation sets in at r > 1 (see Fig. 10). Be- 
cause the system is discrete, there is necessarily a gap 
between the ground state manifold and the low lying ex- 
cited states of the RBIM. In particular. Fig. [Toj suggests 
that this gap is smaller on odd lattices than on even lat- 
tices, so that at a given 'temperature' r there a higher 



chance of that an odd lattice will be in an excited state. 
The origin of this discrepancy remains an open problem. 

Given the maximum possible threshold is set by the 
Nishimorit point at ~ 10.9% [15 , the addition of degen- 
eracy in the matching algorithm improves the threshold 
by a modest but statistically significant margin. 



IV. WHAT CAN WE REALLY SAY ABOUT 
THE RBIM PHASE DIAGRAM? 

The most precise estimate for Pco, the T = phase 
transition in the RBIM, is currently the value reported 
in [15 , and this result has been applied in the study of 
the classical 2D RBIM. The observation of the huge de- 
generacy in the ground state manifold of the RBIM raises 
the question of how to ensure that the ground states (i.e. 
minimum- distance matchings) selected by Edmonds' al- 
gorithm are sampled fairly from the manifold. 

As an elementary test, we applied two freely available 
implementations. Blossom IV [23 and Blossom V |26| , 
to the simple matching problem depicted in Fig. [S) per- 
muting vertex labels through all 24 permutations to re- 
order the presentation of the complete graph. Regardless 
of presentation. Blossom IV always returned the least de- 
generate matching, {{1, 2}, {3, 4}}, regardless of reorder- 
ing. On the other hand, depending on the vertex labels. 
Blossom V returned both possible matchings: half of the 
permutations resulted in matching {{1, 2}, {3, 4}}, and 
the other half yielded the alternative, {{1, 3}, {2, 4}}. 
This test serves to demonstrate that different imple- 
mentations of Edmonds' matching sample amongst the 
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ground states differently. 

To compute Pco accurately, we need to sample the 
ground states of the RBIM fairly, i.e. with a probability 
distribution proportional to the path matching degener- 
acy. To see this, we define an observable Oe,e' for the 
classical RBIM to be the indicator function determining 
whether E -\- E' is in the trivial homology class: 



O 



E,E' 



if E -\- E^ is homologically trivial 

1 otherwise 



As in [15] , Pco is calculated from the scaling behaviour of 
Pfaii = {Oe,e'), which is given by: 

E 

= lim P{E) (9£;,^e-^^^'- /Z 

E Vcr 

E aeg.s. 

= ^P(^) Oe,e'/Ve (17) 

E E'eg.s. 
E M' 

where Dm' is the matching path degeneracy of the 
matching M', Ve = ^M'\dE-^'M ground state 

(g.s.) degeneracy of the RBIM instance E^ the grain 
boundary E' is interchangeable with a given Ising spin 
configuration <j (as described above, and in [10 ), and 



He,(j 



(ij) 



(3JijU 



(19) 



Also we have used the fact that Oe,e' = Oe,m' is con- 
stant over all minimum weight matching paths E' con- 
sistent with a given matching M^ Note that in the limit 
1/ ^ oo, the sum over E is made unnecessary by ap- 
pealing to an ergodic argument: (O) = {0)e-, for any E 
whose errors are uniformly distributed. 

The final result, Eq. (18), establishes the distribution 



from which we should sample the matchings: P{M\E) = 
Dm/'^e- This result suggests that the direct use of exist- 
ing implementations of Edmonds' perfect matching may 
not give accurate estimates of Pco . 

In principle, for a finite lattice, one could enumerate 
all matchings, compute the matching-path degeneracy 
for each, then evaluate Eq. (18) explicitly, however this 



is computationally inefficient. We do not know of a com- 
putationally efficient way to ensure a fair sampling. 



We note however, that within the ground state mani- 
fold (i.e. taking < r < 0.5), the error correction thresh- 
old is bounded by about 0.105 (see Fig.[9|. This suggests 
that Pco < 0.105. That this is statistically less than the 
Nishimori probability Pc gives support to the notion that 
the phase diagram is indeed reentrant, as has been dis- 
cussed elsewhere. However we propose that an accurate 
value of Pco may be higher than that reported in ^15^ . 
V. CONCLUSION 



We have shown using both numerical studies and an- 
alytic results that the simultaneous expression of both 
loss and computational errors can be dealt with in the 
context of topological memories. The trade-off between 
these processes is rather gentle. At low loss rates, our 
analytical results indicate that, because superplaquettes 
embody more physical qubits, the principle effect of losses 
is a moderate, linear increase in the effective error rate 
per plaquette. At higher loss rates, changes in the struc- 
ture and connectivity of the superplaquete lattice become 
significant. 

We have also shown that the classical algorithms used 
to implement topological quantum error correction can 
be improved by accounting for degeneracy in matching 
algorithms. This can take various forms; either in mod- 
ifying edge weights between neighbouring superplaque- 
ttes on an irregular lattice as discussed in pT., or by 
counting the number of degenerate shortest-paths be- 
tween syndromes, as discussed here and in [17]. Both 
approaches may be implemented with little overhead to 
existing matching routines. 

This work also highlights the attention that needs 
to be given to degeneracy when calculating thresholds 
in the RBIM and any related classical statistical me- 
chanics simulations relying on implementations of Ed- 
monds' matching algorithm. To accurately predict Pco, 
minimum-distance matchings should be sampled fairly 
from amongst the hugely degenerate ground state man- 
ifold. To our knowledge, achieving this remains an open 
problem. We hope that this avenue of research will shed 
new light on the RBIM. 
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